Event-driven Monte Carlo algorithm for general potentials 
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We extend the event-chain Monte Carlo algorithm from hard-sphere interactions to the micro- 
canonical ensemble (constant potential energy) for general potentials. This event-driven Monte Carlo 
algorithm is non-local, rejection-free, and allows for the breaking of detailed balance. The algorithm 
uses a discretized potential, but its running speed is asymptotically independent of the discretization. 
We implement the algorithm for the cut-off linear potential, and discuss its possible implementation 
directly in the continuum limit. 
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The event-chain algorithm [l|, a non-local, rejection- 
free Markov-chain algorithm for hard-sphere systems, has 
proved considerably faster than the local Monte Carlo 
algorithm. It has allowed us to show that two-dimensional 
melting in hard disks proceeds via a first-order liquid- 
hexatic transition 0|. 

In the present paper, we generalize the event-chain 
algorithm to the case of general potentials, greatly ex- 
tending the scope of the original method. As for event- 
driven molecular dynamics (MD), we discretize the po- 
tential (here on an energy scale A^;), but unlike the MD 
algorithm [3, 01 , simulations at arbitrary small A e are 
feasible. 
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Figure 1: Event-chain algorithm for hard disks ri)|. Left : 
Collective move, with a total displacement ^ as a fixed para- 
meter of the algorithm. The move is never rejected. Right : 
Return move, demonstrating microreversibility. 



In the event-chain algorithm a randomly chosen 
disk is moved in a straight line by a distance I, if this 
generates no overlap for all £' < £, or else until it hits 
another particle i'. The latter particle, i' , then moves for 
the remainder of the distance £ or, again, until it hits 
yet another particle, etc. (see Fig. [ij. The event-chain 
move is micro-reversible, and it satisfies detailed balance 
if all the displacements are in the same direction. For 



hard disks, it is about two orders of magnitude faster 
than the local Monte Carlo algorithm Q , and it compares 
favorably with the MD methodQ. 

In the event-driven Monte Carlo algorithm, we consi- 
der a discretized potential V{r), with discontinuities such 
that V{r) and thus the total potential energy are mul- 
tiples of a given energy step A^ . From an initial confi- 
guration of energy E, a randomly chosen particle i is 
moved in a straight line until the total displacement, at 
energy E only, equals a fixed value i or else until it "hits" 
another particle i' . This means that the energy of the sys- 
tem would exceed E if the particle i moved any farther. 
The latter particle, i' , then moves for the remainder of 
the distance £ or, again, until it hits yet another particle, 
etc. (see Fig. [2]). In the absence of boundary effects, the 
move is guaranteed to terminate. It is again rejection- 
free and micro-reversible. Detailed balance may again be 
broken and the algorithm can be run with moves, say, in 
the +x and +y direction only. 
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Figure 2: Event-driven Monte Carlo displacement for N = A. 
Particle i moves from its initial configuration until its collision 
event with particle i' (potential of Eq. ([T]) with A_b ~ 1/4 and 
E — 3/4). Only displacements at energy E are discounted 
from £. 

The event-driven Monte Carlo algorithm can be im- 
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plemented by computing the intersection points of the 
potential discontinuities with the trajectory of the par- 
ticle i. One recovers E{x) by adding up the signs in the 
sorted list of the intersections. This is illustrated in Fig. [2] 
for the potential defined in Eq. In Fig. [31 we consider 
particles interacting with the potential 



Vir) = 



if r > 1 

1 — AeIt/Ae] otherwise 
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The root -Econt (a;root ) = E oi the continuous potential is 
uniquely determined via a decision tree in subsequent 
C°° intervals of E{x). For finite A^, one must only iden- 
tify the discontinuity steps of the potential in the interval 
Xroot ±const[(?i?/9a:(a;root)]^^- We have implemented the 
algorithm for very large values of l/A^ (with I oc A^), 
and achieved constant scaling of the algorithm. 

Although the event-driven Monte Carlo algorithm can 
be implemented for arbitrarily small A^;, and does not 



for A^; = 1/10. We compare the pair-correlation function 
of the event-driven Monte Carlo algorithm for a large sys- 
tem with N = 128^ particles to the results obtained with 
a local Monte Carlo algorithm in the same system. With 
the local Monte Carlo algorithm, only configurations with 
potential energy E are used. There is no doubt about 
the correctness of the algorithm and the implementation. 
With this simple implementation, we reached 5 x 10^ 
collisions per hour on a 3 GHz workstation. 
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Figure 3: Pair-correlation function computed with the lo- 
cal and the event-driven Monte Carlo algorithms (A'^ = 128^, 
box of size 128^ E = iVAfi, Ab = 0.1). Data agree, which 
confirms the correctness of the algorithm. 

In the naive implementation suggested in Fig. [2l the 
number of steps between events scales oc 1/A^;, and the 
algorithm becomes very slow in the limit A^; — > 0. Unlike 
for event-driven molecular dynamics Q, this difficulty 
can be overcome, and the speed of the algorithm remains 
constant for A^; — >■ 0. The algorithm is straightforward 
for a convex potential with finite support, as for example 
the continuum limit of Eq. ([T]) , given by 



if r > 1 

1 — 7' otherwise ' 
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for which E{x) is continuous and piecewise C°° while 
dE/dx is piecewise monotonously decreasing (see Fig.|3]). 



Figure 4: Energy E{x) for the configuration of Fig.[51 for dif- 
ferent discretizations, and in the continuum limit, and deriva- 
tive dEcont/dx which is piecewise monotonously decreasing. 



slow down then, it is interesting to consider the Ae — 
limit. A good model for this is to suppose that rather 
than being fixed, A^^ is sampled, between events, from 
a distribution with very large mean. The collision rules in 
Fig. [2] then become probabilistic. They are readily analy- 
zed in simple cases. The event itself, which takes place at 
energy E, is initiated and ended by particles coming in 
from and moving out to lower energies. This event may 
involve a sizable number of particles and thus become 
quite complex. It is unclear whether truly general sto- 
chastic rules for the collisions and the use of £ can be 
obtained. 

In conclusion, in this paper, we have generalized the 
successful event-chain algorithm from hard spheres to the 
case of general potentials and provided details on its im- 
plementation. This algorithm is an alternative to the lo- 
cal Monte Carlo algorithm Q and, even with a simple 
implementation, should outperform it in most applica- 
tions. Moreover, arbitrarily small discretization steps can 
be handled efficiently, which is not the case in the event- 
driven MD ^MM- 
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